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;_( I Abstract: Fresnel diffraction calculation on an arbitrary shape surface is 

Cu . proposed. This method is capable of calculating Fresnel diffraction from 

a source surface with an arbitrary shape to a planar destination surface. 
Although such calculation can be readily calculated by the direct integral 
of a diffraction calculation, the calculation cost is proportional to 0{N^) in 
one dimensional or 0{N'^) in two dimensional cases, where A^ is the number 
of sampling points. However, the calculation cost of the proposed method 
Q ■ is 0{NlogN) in one dimensional or 0{N^logN) in two dimensional cases 

using non-uniform fast Fourier transform. 



in 



OCIS codes: (0070) Fourier optics and signal processing; (090.1760) Computer holography; 
(090.2870) Holographic display; (090.5694) Real-time holography; (090.1995) Digital holog- 
C^ 1 raphy. 
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1. Introduction 

Diffraction calculations such as Huygens diffraction, Fresnel diffraction and angular spectrum 
method, are important tools in wide-ranging optics |IT]|2], ultrasonic 13], X-ray iHJ and so forth. 
Its applications in optics include computer-generated-hologram (CGH) and digital holography 
Is), phase retrieval, image encryption and decryption and so forth. 

Fast Fourier transform (FFT)-based diffraction calculations according to their convolution or 
Fourier transform form are used in these applications. The FFT-based diffraction calculations, 
however, can only be applied to planar surfaces in parallel. In order to apply the methods to 
a non-parallel planar surface, many methods have been proposed: for example, non-parallel 
Fresnel diffractions JS] |7] |8] |9] and non-parallel angular spectrum methods ifTol [TTl [121 [131 . 
Unfortunately, these non-parallel diffractions are limited to planar surfaces. If we calculate a 
source surface with arbitrary shape using non-parallel diffractions, we need to approximate 
the arbitrary shape surface with many small non-parallel planar surfaces, and then, we need to 
calculate the non-parallel diffraction per the small non-parallel planar surfaces. 

In this paper, we propose Fresnel diffraction calculation on an arbitrary shape surface. With- 
out approximating an arbitrary shape surface with small non-parallel planar surfaces, this 
method is capable of calculating Fresnel diffraction from a source arbitrary shape surface to 
a planar destination surface. Although such calculation can be readily calculated by the direct 
integral of a diffraction calculation, the calculation cost is proportional to 0{N^) in one dimen- 
sional or 0{N^) in two dimensional cases, where A^ is the number of sampling points. However, 
the calculation cost of the proposed method is 0{N\ogN) in one dimensional or 0{N^\ogN) 
in two dimensional cases using non-uniform fast Fourier transform. 

In Section 2, we describe the arbitrary shape surface Fresnel diffraction. In Section 3, we 
present the numerical results. Section 4 concludes this work. 

2. Arbitrary shape surface Fresnel diffraction 

Let us begin with Huygens diffraction. Huygens diffraction fl] on a planar surface is ex- 
pressed as: 

"2(X2) = TJ- / /m/(xi)mi(xi) 2 "^^1' <^1) 

where, mi (xi) and M2(x2) are planar source and destination surfaces, xi and X2 are the position 
vectors on the source and destination surfaces, X and k are the wavelength and wave number 

of light, and r — J \x2 — xi\^ + Zq, where zo is the distance between the source and destination 
surfaces. m/(xi) is the incident wave to the source surface. If the incident wave is used as a 
planar wave that is perpendicular to the optical axis, we can treat as m/(xi) = 1. 

We expand the Huygens diffraction to a source surface with arbitrary shape. As shown in 
Fig. [T] the source surface with arbitrary shape Mi(xi,iii) is defined by the displacement di ~ 
d\ (xi) at the position xi. Note that when the arbitrary shape surface is uniform-sampled, the 
corresponding coordinate xi is non-uniform-sampled depending on the slope of u\(ii.\,d\) at 
the position xi. Huygens diffraction on an arbitrary shape surface is expressed by: 

m(X2) = -^ j I M/(xi,rfl)Mi(xi,rfi) 2 ^Xi, (2) 



Here, if the incident wave M/(xi,iii) is used as a planar wave that is perpendicular to the 
optical axis, the incident wave is expressed as ui{x\,d\) = t\p{ikd\). 
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Fig. 1. Diffraction calculation between source surface with arbitrary and planar destination 
surface. 



Applying Fresnel approximation to Eq.Q using ro ^ zo~ di, we can obtain the following 
approximation: 
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(4) 



And, we also approximate r^ ~ Zq in the integral of Eq.(|2]i. Therefore, we obtain Fresnel 
diffraction on an arbitray surface: 



«(X2) = TT — / / M/(Xl,lil)Ml(Xl,rfl) 

2 2 

exp{ik{-di + -J-))exp(-!l )exp{-^)dxi 

Zro ro Zro 



(5) 



The above Eqs.(|2]i and (|5]l can be treated an arbitrary shape surface. We can readily calculate 
by the direct integral with regard to these equations; however, the calculation cost is 0{N^) in 
one dimensional or 0{N'^) in two dimensional cases, where A^ is the number of sampling points, 
because we cannot calculate them using Fourier transform. 

In order to obtain the Fourier form of Eq.(|5), we approximate the third exponential term in 
the integration as follows: 



2 2 

exp(il-^)«exp(/fc-^). 
2ro Zzo 



(6) 



Eventually, we obtain the following equation: 

SXp{ik{zo + 2^)) 
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(7) 



exp[ik )dxi. 

zo~di 



Because the coordinate xi is sampled by the non-uniform sampling rates, instead of (uni- 
form) Fourier transform, we can calculate the above equation using non-uniform Fourier trans- 



(NUFT): 








exp{ik{zo + ^)) 
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(8) 



where, NUF [•] denotes NUFT. NUFT of a function /(xi ) is defined as: 



F(X2)=NUF 



/(xi 



/ / /(xi)exp(-;7rxiX2)c/xi 



(9) 



Although the form of NUFT is similar to that of uniform Fourier transform, the coordinate 
xi is non-uniform- sampled and X2 is uniform-sampled, unlike uniform Fourier transform. For 
the numerical implementation of Eq.(|9]l, it is necessary to use non-uniform fast Fourier trans- 
form (NUFFT) which has the complexity of 0{NlogN). Many methods for NUFFT have been 
proposed over the course of the past twenty years or so [ll4] [TS] [T6l (TT) . NUFFTs are based on 
the combination of an interpolation and the uniform FFT. In this paper, we used L. Greengard 
and J. Y. Lee's NUFFT |.17J. For more details, see ifTTI . 

3. Result 
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Fig. 2. Intensity profiles of diffraction results from a source surface composed of four small 
planar surfaces with 128 points, which are tilted —30°, —50°, 0° and +50°to xi, respec- 
tively, (a) source surface (b) diffraction result by Eq.(|2] (c) diffraction result by Eq.(|5) (d) 
diffraction result by Eq.(|8) (e) Absolute error between (b) and (d). 



Let us examine our method using two source surfaces with arbitrary shapes in one dimension: 
a surface composed of four tilted planar surfaces, and quadratic curve. 



We used the distance zo = 1 m, the wavelength of 633 nm, and the number of sampling points 
on source and destination N = 1, 024. The sampling rates on the source and destination surface 
are p = 10/i m. We used a planar wave as the incident wave that is perpendicular to the optical 
axis. 

Figure |2] shows a source surface composed of four small planar surfaces with 128 points, 
which are tilted —30°, —50°, 0° and +50° to xi, respectively. The horizontal axis in (a) in- 
dicates the position on xi in metric units. The horizontal axes in (b)-(e) indicates the position 
on X2 in metric units. The sampling rates on these small planar surfaces are p, however, the 
sampling rates on xi are |/9Cos(— 30°)|, |/?cos(— 50°)|, |/?cos(0°)| and |/?cos(+50°)|, respec- 
tively. The destination planar surface is not inclined to X2. Figures |2](b)-(d) show the intensity 
profiles of the diffraction results by Eq.(|2), Eq.(|5) with direct integral and our method (Eq.(|8)), 
respectively. Figure|2](e) depicts the absolute error between (Eq.©) and Eq.®. The absolute 
error falls into within approximately 0.025. 
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Fig. 3. Intensity profiles of diffraction results from quadratic curve surface, (a) source sur- 
face (b) diffraction result by Eq.(|2) (c) diffraction result by Eq.l|5j (d) diffraction result by 
Eq.© (e) absolute error between (b) and (d). 

Figure |3] shows a source surface with quadratic curve. The sampling rate on Xi according to 
the source surface depends on the slope of the quadratic curve. Figure|3](b)-(d) shows the inten- 
sity profiles of the diffraction results by Eq.Q, Eq.(|5) and our method (Eq.®), respectively. 
Figure |3](e) depicts the absolute error between Eq.© and Eq.®. The absolute error falls into 
within approximately 0.005. The primary factor of these absolute errors in Figsj2]and|3]is the 
approximations by Eqs.© and (|6]l. 



4. Conclusion 

We proposed Fresnel diffraction calculation from a source surface with arbitrary shape to a 
planar destination surface, using one NUFFT. We have had only a direct integral for calculating 



such diffraction thus far. Unfortunately, it is very time consuming because the calculation cost 
takes 0{N-^) in one dimensional or 0{N'^) in two dimensional cases. In contrast, the calculation 
cost of our method is 0{NlogN) in one dimensional or 0{N-- \ogN) in two dimensional cases 
using NUFFT. The method is very useful for calculating a CGH from a three-dimensional object 
composed of multiple polygons or arbitrary shape surfaces. In our next work, we will show the 
fast calculation of a CGH from such 3D objects using this method. 
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